1 Introduction 

In the present paper, we consider an interface between two liquid layers subject to hori- 
zontal harmonic vibrations, of amplitude a and frequency uj. 

In experimental works by Bezdeneznykh et al. |]T| and by Wolf for a long horizontal 
Qi^'. reservoir filled with two immiscible viscous fluids, an interesting phenomenon was found 
^ \ at the interface: the horizontal vibrations lead to the formation of a steady relief. This 
formation mechanism has a threshold nature; it is noteworthy that such a wavy relief appears 
on the interface only if the densities of the two fluids are close enough, i.e. it does not appear 



for the liquid/gas interface (free surface). The interface is absolutely unstable if the heaviest 
fluid occupies the upper layer; i.e., the horizontal vibration does not prevent the evolution 
of Rayleigh- Taylor instability, in contrast to the vertical one which under certain conditions 



< 



^ ! suppresses its evolution. A theoretical description of this phenomenon was provided by 

^ \ Lyubimov & Cherepanov within the framework of a high frequency (of the vibration) 
O ' 

■ approximation and an averaging procedure; they found that a horizontal vibration leads 
>>• 

(-| . to a quasi-equilibrium state i.e., a state where the mean motion is absent but the interface 

oscillates with a small amplitude (of the order of magnitude of the cavity displacement) with 

^ _ respect to the steady relief. In particular, they found that the rehef (in the case of a heavy 

^ ■ fluid occupying the bottom layer) with finite wavelength is not possible for any thickness of 
O ! 

00 ■ fluid layers. In fact, the relief with finite wavelength arises only for considerably thick layers 

O ■ 1/2 

On ' of height h > [3a/ {g (pi — P2))] , where a is the coefficient of surface tension, pi and p2 
0\ ' 

■ are the densities of lower and upper layer, respectively; g is the gravity acceleration. In the 

o : 

. approach [Q, two parameters were assumed to be asymptotically small simultaneously: (i) 



^ ' the dimensionless thickness of the viscous skin-layers 5 = h ^Ju/to, v being the kinematic 

Ph. ^ 

viscosity and (ii) the dimensionless amplitude of the vibration e = a/h. But in the limiting 

^ ■ case e — s> 0, the possibility of description of parametric resonant effects is absent and only 

H ; 

5t , the basic instability mode (Kelvin-Helmholtz, of two counter flows) remains. 

In this investigation, the condition e — > is dropped and the vibration frequency is 
considered as relatively low, thus the averaging method could no longer be applied. 

2 Problem formulation 

Let us consider the system of two immiscible, incompressible liquids filling rectangular 
cavity of length L and height h. In the state of rest the heavy liquid (of density pi) occupies 



the bottom region of height hi, and the light Hquid (of density P2) - the upper region of 
height /i2 {h = hi + h2)- We choose Cartesian coordinate system in such a way that the 
X, y-axis he in horizontal plane, the 2;-axis is directed vertically, z — Q corresponds to the 
unperturbed interface. Let the cavity perform harmonic oscillation along the x-axis, with 

— * 

the amplitude a and frequency a;, and let 7 be the unit vector along the z-axis, while j along 
the vibration axis (Fig. 1). 

The waves are generated on the interface due to such oscillation, and propagate from 
the vertical end walls to the central part of the cavity. But, if viscosities of the liquids are 
not too small, such a waves are damped on relatively short distances from the endwalls. 
For example, for transformer oil and glycerine in Earth gravity field and angular frequencies 
of the vibration of the order of magnitude lO^Hz, the damping length Id is of the order 
of magnitude of 1 mm. In the following, we assume that the condition L 3> is satisfied 
and we consider only the central part of the cavity, far away from the endwalls, which the 
interface waves generated near the endwalls do not reach. Therefore in the basic state the 
interface could be considered as plane and horizontal. 

However, the liquid in the central part of the cavity can't stay at rest, and the velocities 
of the liquid motion in the bottom and upper layers must be different. Indeed, since the 
densities of liquids are different, then pressures gradients would be different if the velocities in 
the layers are equal. The latter is in contradiction with the normal stress balance condition 
at the interface. The only exception might be the case of liquids which stay in rest in 
laboratory reference frame, but such a situation is impossible because of the incompressibility 
conditions. To prove this, let us consider the cross-section of the cavity by the vertical plane 
— + S+, attached to the cavity (Fig. 1). Then hquid volume V — Vj— + K+ to the 
right (for example) of the cross-section is constant. By virtue of the incompressibility the 
latter means that in the cavity frame the full liquid consumption through the must be 
zero: 



where ^ is the vertical coordinate of the interface. We stress that the FULL consumption 
across any vertical cross-section is zero. At the same time the consumptions of individual 
liquids can be non-zero if the flow is unsteady because the interface is not at rest near the 
vertical endwalls, therefore the volume of each liquid to the right of the cross-section is 
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Fig. 1. 



time-dependent. If the liquids stay in rest in the laboratory frame, then in (2.0) the velocities 
and V2x are equal and (2.0) could not be satisfied. Thus, despite of the statement "in long 
cavity the interface waves generated near the endwalls do not penetrate into the central part 
of the cavity" , the motion in the central part is always present and all of the above discussion 
gives the opportunity for the approximation of the infinite horizontal layer. Moreover, the 
condition (2.0) must be, of course, juncted. 

The concrete flow in the inviscid case is obtained in Section 3, the flow in viscous case 
- in Section 4.2. We stress as early as in this Section that the basic flow, as is obvious 
from (2.0), is of counter-flow type, that makes us to expect the rise of Kelvin-Helmholtz 
instability, while the periodic character of the flow allows to expect the parametric resonance 
effects. 

In the cavity frame, the Navier-Stokes equations for the motion could be obtained from 
the standard equations written in the laboratory frame by considering the total acceleration, 
i.e. 

~91 ~^ ~91 + o.'^^j cos ujt. 
And so, the full Navier-Stokes equations read 



^ + (v/sV)vp = — -"^Pp + vp^vp -g^ + au'^j cos ut 
at pp 



(2.1) 



and the continuity equations for the two layers read 

div vp — 



(2.2) 
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where the subscript = 1, 2 refers to the lower and upper layers, respectively; the other 
notation being conventional. 

No-slip boundary conditions for the velocity are imposed on the rigid walls of the cavity 



at z^-hi: vi = (2.3) 

at z = /i2 : = (2.4) 

and on the interface z — ^{x, y, t) , the following conditions are satisfied: 

- stress balance 

(Pi - P2)ni = {all^ - a'^k)nk + arii div ft (2.5) 

- velocity continuity 

vi = V2 (2.6) 



kinematic condition 



+ {v,V)^ -v,-^ (2.7) 



In (2.5 - 2.7), a is the coefficient of the surface tension, are the viscous stress tensors: 



~ dxk dxi 



The integral condition of zero full consumption reads 



j vi-jdz + J V2-jdz^0 (2.8) 



Vx : 

-hi e 



3 Inviscid approximation 

After omitting the viscous term the Navier-Stokes equations (2.1) read: 

+ ivfjV)vf3 = - gj + auj^ J cos ujt (3.1) 

at pp 

The continuity equations remain unaltered: 



div (3.2) 
4 



The impermeability boundary conditions are imposed on tlie rigid walls instead of the no-slip 
conditions: 

at z — —hi : Vi ■ J — (3.3) 
ai z — h2 : ^2 • 7 = (3.4) 

On the interface z = ^{x,y,t) the following conditions are imposed: 

- normal stress balance 

Pi — P2 — CK div n (3.5) 

- continuity of normal components of velocity 

Vi ■ n — V2 ■ n (3-6) 

- kinematic condition 

§ + {viV)i = iTi • 7 (3.7) 

The zero full consumption condition (2.8) remains unaltered. 
The problem (3.1) - (3.7), (2.8) admits a solution of the form: 

Vp = Upjsmijjt (3.8) 

where 

Ui = au)- — ^ ; — - — aUi and U2 = —auj- — ^ ; — - — aU2 (3.9) 

hiP2 + h2Pl hiP2 + /l2Pl 

The corresponding pressure field is given by the expression 

P3 — —pggz + auj P1P2-, ; xcosujt (3.10) 

hip2 + h2pi 

The above solution corresponds to a plane-parallel, unsteady, counter flow which keeps the 
interface fiat ( ^ = ). Note that in the case of equal densities, Vi and V2 tend to zero, 
i.e. fluids stay at rest in the cavity frame. In the case P2 <S Pi , the lower fluid stays at rest 
in the laboratory frame {Ui — auj ). 

Thus, in the inviscid case we deal with the stability problem for the interface between two 
counter flows, and the difference in the flows velocities is a periodic function of time. Linear 



stability problem could be reduced to the ordinary differential equation for the amplitude 
^{t) of the interface displacement from quasi-equilibrium horizontal position: 



(Fi + Fs)^ + 2ik'^{FiUi + F2U2) sincut+ 

+ (pi - p2)gk + i{FiUi + F2U2) ku cos tot - P{FiUl + Faf/I) sin^cut) = 0, (3.11) 

with Fi = pi coth(/c/ii) and F2 = P2 coth.{kh2) 

In (3.11), k is the wavenumber of the disturbances, which determines the periodicity of the 
solution along vibration axis. It is straightforward to show that stability problem admits the 
analogue of Squire's theorem The latter gives the possibility of considering only plane 
(2D) disturbances. 

It is convenient to eliminate from (3.11) the term which contains the first order derivative 
in ^ . To do so, we make the change of variable 

^(t) = y(t)e'*W (3.12) 

where 

k F,U, + F2U2 . 
<I> = — coscut 3.13 

iO F1+F2 ^ ^ 

(since $ is real, then ^ and Y are equal via a modulus, i.e. ^ and Y are equivalent from the 
point of view of stability). This results in the standard Mathieu equation for Y: 

^ + iA-QcosH)Y = (3.14) 

The equation (3.14) is in dimensionless form, with the following reference quantities for time 
and length: 

f = 1/iu and l* = k = H {g{pi - P2))Y'^ 

where Ic is the capillary length, and with the following notation: 

_ AB,e pcothjkHi) coth{kH2) {Hi + ff2)'(p - 1)' . . 

^ Wei {pcoth{kHi) + coth{kH2)y {H1 + H2PY ^' ^ 

^ k{l + e) p-1 

Wei pcoth{kHi) + coth{kH2) ^' ' 
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In formulas (3.15) and (3.16), Wei = ^"^^c/ 9 is Weber number based on the capillary length, 
p = pi/p2 is densities ratio. Hi = hi/l^. and H2 = h2/lc are dimensionless layer thicknesses. 
We also introduced the dimensionless parameter B^, characterizing the vibrations intensity: 



By 



a^b? ( pi - -^/^ 
ga 



The solutions of the equation (3.14), which correspond to the neutral boundary, form two 
classes: having period 2-7? (harmonic disturbances) and YL, having period vr (subharmonic 
disturbances). The solutions of subharmonic type are antiperiodic, that means that they 
change the sign for a shift of tt. 

We plot the boundaries of the instability regions (obtained by numerical integration of 
the Mathieu equation (3.14)) in terms of the two parameters, B^ and k. The neutral curves 
shown in Fig. 2, 3 correspond to Weber numbers equal respectively to Wei = 10 and 
Wei = 100 {Hi = H2 = l,p = 2). 

The left curve in this Figures bounds Kelvin-Helmholtz instability region. In the referred 
case Hi = H2 = 1 the most dangerous are longwave disturbances with k = 0. In was 
investigated the case Wei —>■ 00 with finite heights of the layers. It was found that if 
Hi = H2 = H then the transition from the longwave instability to the finite wavelength 
instability occurs aX H = H^ = y^. It is straightforward to analyse the case of finite Wei. 
Note, that if k is small, then Q and A are small also, and could be expanded into power 
series in k. The result is that ai H < H^ the longwave instability takes place at any values 
of Wei, but at if > if * only for 

In the instability region discussed, the instability mechanism is, in fact, the same as in 
the standard case of Kelvin-Helmholtz instability: the velocity grows up and the pressure 
falls down over the elevations of the interface. The consequence of this is the reducing of the 
effective stiffness of the system and, starting from some threshold, the rise of the instability. 
As shown in , in the case of high frequencies of the vibration the nonlinear developpment 
of such instability results in the formation of a quasi-stationary wave relief on the interface. 

However, periodic changes of the velocity can result not only in the average effect, but in 
the resonance amplification of eigenoscillations. Inasmuch as in the absence of the dissipation 



the eigenoscillations arc not damped, then the instabihty will take place at infinitely small 
amplitude of the vibration (if the synchronism condition is satisfied). On Fig. 2, 3 the regions 
of parametric instability approach the abscissa as narrow ends ("tongues"); the points kn 
at which the n-th region of instability is in contact with the /c-axis can be evaluated 
from the equation A — n^,n — 1,2,.... As the vibration frequency grows, i.e. as Weber 
number increases, the points A;„ shift to the shortwave region (as seen from Fig. 2,3) and 
in the limit Wei oo they are displaced into infinity. Thus, at high frequencies of the 
vibration the parametric instability can take place only for shortwave disturbances but, for 
such disturbances, the viscosity could not be ignored. It has to be expected that for big 
values of Weber number only the quasi-static instability must be observed, with mechanism 
of excitation which is weakly sensitive to the viscosity. 

4 2D linear stability problem for viscous fluids 

To simplify the following presentation, we assume in this Section that the two layers are 
of equal height hi — h2 — h. In addition, we non-dimensionalize the problem using the 
scales 

t* — l/cu; I* — h; u* — acu; p* — p2hajuP 

4.1 Problem formulation 

The following dimensionless equations replace the equations (2.1) - (2.8): 

^ + A{vpV)vp = -RpVpp + n-p^Avp - GoA-^^ + kcost (4.1) 

div vp^O (4.2) 

at ^ = -1 : vi=0 (4.3) 

at ^ = 1 : V2^0 (4.4) 

at 2; = ^{x, y, t) : \p]ni = {^^^pa^^^ — Q:2^'^ik)^k + A~^We2^ni div n (4.5) 
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[^;] = 



(4.6) 



^^ + (^iV)e = ^^i-7 (4.7) 
i 1 

Vx : J vi-j dz + j V2-jdz^0 (4.8) 

Here and below, the quantity jump across the interface is denoted by square brackets, for 
example \p\ — Pi — P2, and the following parameters are introduced: 

A — ah"^; ilp — h^uu^^; Go — g/ {huj'^)] We2 — p2h^ou^/a; p — p\/p2\ Ri — l/p and R2 — 1. 

So, the system of equations governing the fluid motion contains 7 dimensionless parameters: 
A is dimensionless amplitude, VLjS are dimensionless frequencies, the parameter Go, We2 is a 
new Weber number, k is the wave number and p is the densities ratio. 

4.2 Solution method 

The system of equations (4.1) - (4.8) admits a solution of the form 

Vp = Up{z)]ei^v{it)+C.C. (4.9) 

Pp = -GoA-^R^^z + x{SeyLY>{it) + C.C.), S = const. (4.10) 

The above solution corresponds to the plane-parallel unsteady counter flow which keeps the 
interface flat ( C = ) (to be compared with (3.8)- (3. 10) in the inviscid case). The function 
Ui3{z) could be found analytically, but its expression is too complicated, so it was determined 
numerically. Figure 4 represents the velocity proflles at four time instants: (a) t = 0, (b) 
t = 7r/2, (c) t = 71, (d) t = 37r/2. The direction of the flow changes every half period of the 
external forcing. Note the thin boundary layers at both sides of the interface and near the 
walls ( Figs. 4a, 4c). 

For the reasons discussed in Section 3, we investigate the stabihty of the basic state 
(4.9), (4.10) with respect to periodic (in time and space) 2D disturbances. If the vibration 
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frequency is high enough, then thin boundary layers occur at the interface (on both sides) and 
near the two horizontal walls, making the task of numerical solution particularly complicated. 

After linearizing near the solution (4.9) and (4.10), the Navier-Stokes equations for small 
disturbances and pp are: 

^ + A{Vf3V)uf3 + A{ufsV)Vfs = -RpVpp + Q^'Auf, (4.11) 



div up^O (4.12) 
Boundary conditions for the system (4.11) - (4.12) read: 

at z = -1 : ui^O (4.13) 

at z = l: U2^0 (4.14) 

at ^ = : ((1 - p)GoC + We^'A^) + = {^-'pag^ - n^'a^l^)nk (4.15) 



dV 
dz 



(4.16) 



^§ + (W)e = ^^i. (4.17) 

This is valid as long as the deformation ^ is small compared to the wavelength of the 
instability. 

It is easy to see that the problem is uniform in x-direction. This allows us to consider 
only normal disturbances. The linear stability problem is now summarized as: 

+ ikAVpUp + wpA-^ = -ikRppp + - ^')^^ (4-18) 

^ + ikAV,w, = -i?,^ + n-( ^ - e)w, (4.19) 
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ikU3 + ^r— = 

oz 



(4.20) 



at z — —1 : Ui = Wi = 



(4.21) 



at 2; = 1 : U2 — W2 — 



(4.22) 



atz = 0: A-^^ ((1 - p)Go - We^^k^) + b] + 2 (^2 



^^"^-l^rV^Uo (4.23) 



dz 



dz 



^2 



p ikwi + — h C 



dui ^d^Vi\ , (9^2 



9^2 



= 



(4.24) 



^2; 



(4.25) 



H = 



(4.26) 



(4.27) 



All fields could be expanded in time in Fourier series of the form 

00 

t) = M^) exp(mt) + C.C., 

n=— 00 

where n might be taken integer, as well as half-integer. By keeping only a finite number of 
terms of the expansion in both fluids we get 2-point boundary value problems (EVP's, cou- 
pled through the conditions (4.23) - (4.27)), for the system of ordinary differential equations 
for the complex amplitudes pn/s^z), Unpiz), Wnpiz) (note that Cn{z) may be eliminated from 
(4.23)-(4.27)). The coefficients in the equations depend on z, therefore the system could be 
resolved only numerically. 

One important feature of the resulting ODE's systems is that integer and half-integer 
harmonics are not coupled (since the basic flow (4.9) introduce only a shift of ±1 for n 
in (4.18)-(4.27)). So, finally, the corresponding EVP's can be solved independently. The 
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solutions of the EVP's involving integer harmonics have the same period as that of the 
forcing vibration (harmonic, or synchronous case) and the solutions of the EVP's involving 
half-integer harmonics have a period twice as that of the forcing vibration (subharmonic, or 
asynchronous case). 

EVP's, in turn, allow the reduction to initial values problems (IVP's). We focus our 
attention on the synchronous case, since asynchronous solutions were not detected in the 
numerical investigation. The numerical procedure is described in 

5 Numerical results 

The linear stability of the fiat interface has been calculated using 21 basis functions 
fn{z) in both cases (synchronous and asynchronous). Ten to twelve basis functions are 
enough to get sufficiently accurate numerical results if the dimensionless frequency is not 
high {Q <^ 250). The number of basis functions has to be increased if > 250. 

In Figure 5, we show the neutral stability curves which divide the {k, y4)-plane into 
regions of stable solutions, and regions (tongues) of unstable (exponentially growing in time) 
solutions. The values of the nondimensional parameters are p=2, 6*0=0.16 and We2=Q.25. 
Only synchronous solutions were detected. This situation is not unique (for the review 
and discussion see 0). Asymptotic "inviscid" curves are plotted with solid lines, while 
the curves obtained numerically for viscous fluids at different values of the dimensionless 
frequency Qi = Q2 = are plotted with dashed lines (note that Qi = Q2 correspond to 

For the high values of Q (corresponding to small viscosities) the numerical curves ap- 
proach the bottom of the asymptotic tongues, which correspond to an instability threshold, 
but lower and lower (Fig. 5). As Q decreases (or the viscosities increase), the instability 
threshold grows (Fig. 6). We state that viscosity has a weak influence on the Kelvin- 
Helmholtz instability, but strongly damps the parametric instability associated with the 
excitation of capillary-gravity waves. In addition, as seen in Figures 5,6, at some critical 
frequency 100 < < 250, the boundary of Kelvin-Helmholtz instability region unites with 
the boundary of the resonant zone, so that at smaller frequencies there exists a single neutral 
curve corresponding to each frequency value. This means that in the case of rather viscous 
fluids the difference between the parametric and non-parametric instabilities vanishes. 
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FIGURES CAPTIONS 



Fig. 1. Problem configuration. 

Fig. 2. Neutral stability diagram (inviscid case), for Wei — 10. 
Unstable tongues touch the A;-axis at points kn, 
defined by the equation (3.28). The bigger is the 
tongue's number, the narrower it is. 
The tongues shift towards shortwave region as the 
vibration frequency increases (next figure). 

Fig. 3. Neutral stabihty diagram (inviscid case), for Wei — 100. 
The comparison of this Figure and Figure 2 reveals 
that the influence of changes in the vibration frequency 
on Kelvin-Helmholtz instability is very small. 

Fig. 4. Velocity profile at four time instants (viscous case): 
(a) t = 0, (b) t = tt/2, (c) t = tt, (d) t = 37r/2. 
Flow direction changes every half period of the 
forcing vibration. 

Fig. 5. Neutral stability k — A diagram (viscous case). 

Parameters are: p=2, 6*0=0.16 and We2=6.25; 
il is large (250, 360), i.e., viscosity is small . 

Fig. 6. Neutral stabihty k — A diagram (viscous case). Parameters are 
the same as in Fig. 5; Q is small (25, 50, 100), i.e., 
viscosity is large. Viscosity makes weak influence on Kelvin-Helmholtz 
type instability but strongly damps the parametric instability. 
In the case of rather viscous fluids, the difference between 
the parametric and non-parametric instabilities vanishes. 
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